###########################
# Table 1 and Table A4
###########################

rm(list=ls())       # clear objects in memory

#sink("~/Dropbox/COVID Peru/08_replication/02_randomization_inference_economy.txt")

library(ri)         # load the RI package
library(foreign)    # package allows R to read Stata datasets

set.seed(123)

d = read.dta("clean_survey_experiment_covid_peru.dta")
d_eco = d[!is.na(d$t_ind_economia),]
d_eco2 = d_eco[!is.na(d_eco$acuerdo_cuarentena),]

table(d_eco$t_ind_economia)
table(d_eco2$t_ind_economia)

#######################################
# Economy treatment Binary outcome
#######################################

Z <-  d_eco$t_ind_economia      
Y <- d_eco$acuerdo_cuarentena_binary

probs <- genprobexact(Z)  
table(probs)

ate <- estate(Y,Z,prob=probs)      # estimate the ATE
ate

perms <- genperms(Z,maxiter=10000)   # set the number of simulated random assignments

Ys <- genouts(Y,Z,ate=0)    # create potential outcomes under the sharp null of no effect for any unit

distout <- gendist(Ys,perms,prob=probs)  # generate the sampling distribution  based on the schedule of potential outcomes implied by the null hypothesis

dispdist(distout,ate)       # display p-values, 95% confidence interval, standard error under the null, and graph the sampling distribution under the null

round(mean(Y[Z==1]),3)
round(mean(Y[Z==0]),3)
ate

#######################################
# Economy treatment Ordinal outcome
#######################################

Z <-  d_eco2$t_ind_economia      
Y <- d_eco2$acuerdo_cuarentena

probs <- genprobexact(Z)  
table(probs)

ate <- estate(Y,Z,prob=probs)      # estimate the ATE
ate

perms <- genperms(Z,maxiter=10000)   # set the number of simulated random assignments

Ys <- genouts(Y,Z,ate=0)    # create potential outcomes under the sharp null of no effect for any unit

distout <- gendist(Ys,perms,prob=probs)  # generate the sampling distribution  based on the schedule of potential outcomes implied by the null hypothesis

dispdist(distout,ate)       # display p-values, 95% confidence interval, standard error under the null, and graph the sampling distribution under the null

round(mean(Y[Z==1]),3)
round(mean(Y[Z==0]),3)
ate

sink()

